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Abstract 

We introduce a mathematical framework where the statistics of spikes trains, produced by neural networks 
evolving under synaptic plasticity, can be analysed. 
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1 Introduction. 

Synaptic plasticity occurs at many levels of organi- 
zation and time scales in the nervous system (TJ. It is 
of course involved in memory and learning mech- 
anisms, but it also alters excitability of brain area 
and regulates behavioural states (e.g. transition be- 
tween sleep and wakeful activity). Therefore, un- 
derstanding the effects of synaptic plasticity on neu- 
rons dynamics is a crucial challenge. However, the 
exact relation between the synaptic properties ("mi- 
croscopic" level) and the effects induced on neurons 
dynamics (meso-or macroscopic level) is still highly 
controversial. 

On experimental grounds, synaptic changes 
can be induced by specific simulations conditions de- 
fined through the firing frequency of pre- and post- 
synaptic neurons J2] |3J, the membrane potential of 
the postsynaptic neuron [@), spike timing J5] |6] |7) 
(see |0 for a review). Different mechanisms have 
been exhibited from the Hebbian's ones [9| to Long 
Term Potentiation (LTP) and Long Term Depression 
(LTD), and more recently to Spike Time Dependent 
Plasticity (STDP) J6] |7) (see OH \TT\ El for a re- 
view). Most often, these simulation are performed in 
isolated neurons in in vitro conditions. Extrapolat- 
ing the action of these mechanisms to in vivo neural 
networks requires both a bottom-up and up-bottom 
approach. 

This issue is tackled, on theoretical grounds, 
by infering "synaptic updates rules" or "learning 
rules" from biological observations lfl3l Q] Q4] (see 
ifTUl [TTl Q~2] for a review) and extrapolating, by the- 
oretical or numerical investigations, what are the ef- 
fects of such synaptic rule on such neural network 
model. This approach relies on the belief that there 
are "canonical neural models" and "canonical plas- 
ticity rules" capturing the most essential features of 
biology. Unfortunately, this results in a plethora of 
canonical "candidates" and a huge number of papers 
and controversies. Especially, STDP deserved a long 
discussion either on its biological interpretation or 
on its practical implementation ifTSI . Also, the dis- 



cussion of what is relevant for neuron coding and 
what has to be measured to conclude on the effects 
of synaptic plasticity is a matter of long debates (rate 
coding, rank coding, spike coding ?) |[T6l[T7l 

In an attempt to clarify and unify the overall 
vision, some researchers have proposed to associate 
learning rules and their dynamical effects to gen- 
eral principles, and especially to"variational" or "op- 
timality" principles, where some functional has to 
maximised or minimised (known examples of vari- 
ational principles in physics are least-action or en- 
tropy maximization). Dayan and Hausser [18] have 
shown that STPD can be viewed as an optimal noise- 
removal filter for certain noise distributions, Rao and 
Sejnowski HUEQ) suggested that STDP may be re- 
lated to optimal prediction (a neuron attemps to pre- 
dict its membrane potential at some time, given the 
past). Bohte and Mozer |2]] proposes that STPD 
minimizes response variability. Chechik [22 1 relates 
STDP to information theory via maximisation of 
mutual information between input and output spike 
trains. In the same spirit, Toyoizumi et al ll23l l24l 
have proposed to associate STDP to an optimality 
principle where transmission of information between 
an ensemble of presynaptic spike trains and the ouput 
spike train of the postsynaptic neuron is optimal, 
given some constraints imposed by biology (such 
as homeostasy and minimisation in the number of 
strong synapses, which is costly in view of continued 
protein synthesis). "Transmission of information" 
is measured by the mutual information for the joint 
probability distribution of the input and output spike 
trains. Therefore, in these "up-bottom" approaches, 
plasticity rules "emerge" from first principles. 

Obviously, the validations of these theories re- 
quires a model of neurons dynamics and a statisti- 
cal model for spike trains. Unfortunately, often only 
isolated neurons are considered, submitted to input 
spike trains with ad hoc statistics (typically, Poisson 
distributed with independent spikes [23, 24]). How- 
ever, adressing the effect of synaptic plasticity in 
neural networks where dynamics is emerging from 
collective effects and where spike statistics are con- 
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strained by this dynamics seems to be of central im- 
portance. 

This issue is subject to two main difficulties. 
On the one hand, one must identity the generic 
dynamical regimes displayed by a neural network 
model for different choices of parameters (includ- 
ing synaptic strength). On the other hand, one must 
analyse the effects of varying synaptic weights when 
applying plasticity rules. This requires to handle 
a complex interwoven evolution where neurons dy- 
namics depends on synapses and synapses evolu- 
tion depends on neuron dynamics. The first aspect 
has been adressed by several authors using mean- 
field approaches (see [25] and references therein), 
"markovian approaches" [26], or dynamical system 
theory (see |27| and references therein). The sec- 
ond aspect has, to the best of our knowledge, been 
investigated theoretically in only a few examples 
with hebbian learning ll28l l29l [30l or discrete time 
Integrate and Fire models with an STDP like rule 

EDE2. 

Following these tracks, we have investigated 
the dynamical effects of a subclass of synaptic 
plasticity rules (including some implementations of 
STDP) in neural networks models where one has 
a full characterization of the generic dynamics ll33l 
[34l . Thus, our aim is not to provide general state- 
ments about synaptic plasticity in biological neural 
networks. We simply want to have a good mathemat- 
ical control of what is going on in specific models, 
with the hope that this analysis should some light on 
what happens (or does not happen) in "real world" 
neural systems. Using the framework of ergodic the- 
ory and thermodynamic formalism, these plasticity 
rules can be formulated as a variational principle for 
a quantity, called the topological pressure, closely re- 
lated to thermodynamic potentials, like free energy 
or Gibbs potential in statistical physics |35|. As a 
main consequence of this formalism the statistics of 
spikes are more likely described by a Gibbs proba- 
bility distributions than by the classical Poisson dis- 
tribution. 

In this communication, we give a brief outline 
of this setting and provide a simple illustration. Fur- 
ther developments will be published in an extended 
paper. 

2 General framework. 

Neuron dynamics. In this paper we consider a 
simple implementation of the leaky Integrate and 



Fire model, where time has been discretized. This 
model has been introduced and analysed by mean- 
field technics in 1321 . A full characterisation of its 
dynamics has been done in [33 1 and most of the con- 
clusions extend to discrete time Integrate and Fire 
models with adaptive conductances [34]. 

Call Vi the membrane potential of neuron i. Fix 
a firing threshold 9 > 0. Then the dynamics is given 
by: 

N 

Vi(t+1) = 7 V,(t) (1 - ZM(t)])+£V tf ZM(t)]4 

3=1 

(1) 

i = 1 . . . N, where the "leak rate" 7 E [0, 1[ and 
Z[x) — x( x > G) where x is the indicatrix function 
namely, Z(x) — 1 whenever x > 9 and Z(x) = 
otherwise. Wij models the synaptic weight from j to 
i. Denote by W the matrix of synaptic weights. If xt 
is some (time independent) external current. Call 
j(ext) jjjg vec t or f i<r xt _ yy i( ea:t ) are control pa- 
rameters for the dynamics. To alleviate the notation 
we write T = (W, l( ext )). Thus T is a point in a 
N 2 + N dimensional space of control parameters. 

It has been proved in ll33l that the dynamics of 
(TQ admits generically finitely many periodic orbits 
[J. However, in some region of the parameters space, 
the period of these orbits can be quite larger than 
any numerically accessible time, leading to a regime 
practically indistinguishable from chaos. 

Spikes dynamics. Provided that synaptic weights 
and currents are bounded, it is easy to show 
that the trajectories of CO stay in a compact set 
■M = [Vmin, V max ] N . For each neuron one can 
decompose the interval T = [V m in , V max ] into 
To Uli with X = [Vmin, 0[, li = [6,V max \. 
If V E 2q the neuron is quiescent, otherwise it 
fires. This splitting induces a partition V of M, 
that we call the "natural partition". The elements 
of V have the following form. Call A = {0, 1} N , 

where 



E A. Then, V 



Vu} = V Ul X V U2 X • • • X V WN . Equivalently, if 
V E Pu>, then all neurons such that uji = 1 are 
firing while neurons such that cj s ; = are quiescent. 
We call therefore uj a "spiking pattern". To each 
initial condition VsMwe associate a "raster plot" 
w = {u>(f)}££ such that V(i) E P W (t), Vt > 0. 
We write V — > Q. Thus, uj is the sequence of spik- 
ing patterns displayed by the neural network when 
prepared with the initial condition V. We denote 
by [w] t the sequence w(s), u>(s + 1), . . . ,u>(i) of 



As a side remark, we would like to point out that this models admits generically a finite Markov partition, which allows to describe the 
evolution of probability distributions via a Markov chain. This gives some support to the approach developped in [ 26 ] though, our analysis also 
shows, in the present example, that, in opposition to what is assumed in [26], (i) the size of the Markov partition (and especially its finitness) 
does not only depend on membrane and synaptic time constants, but on the values of the synaptic weights; (ii) the probability that, exactly j 
neurons in the network fire and exactly N — j neurons does not fire at a given time, does not factorize (see eq. (1) in [26 ]); (iii) the invariant 
measure is not unique. 



spiking patterns. We say that an infinite sequence 
ui = {i*j(t)}^J^ is an admissible raster plot if there 
exists V e M such that V — > Q. We call E r the 
set of admissible raster plots for the parameters T. 
The dynamics on the set of raster plots, induced by 
(Q~|l, is simply the left shift er r (i.e. ui' = a r u> is 
such that u'(t) = ui(t + 1), t > 0). Note that there 
is generically a one to one correspondence between 
the orbits of (Q]i and the raster plots (i.e. raster plots 
provide a symbolic coding for the orbits) ||33l . 

Statistical properties of orbits. We are inter- 
ested in the statistical properties of raster plots 
which are inherited from the statistical properties 
of orbits of ([TJ via the correspondence v[A\ — 
fj, [{V — > ui, lu G A}}, where v is a probability mea- 
sure on E r , A a measurable set (typically a cylin- 
der set) and p, an ergodic measure on M.. v(A) can 
be estimated by the empirical average tt^ 1 



one has 



i J2t=i xi^r^ £ A) where x is the indicatrix func- 
tion. We are seeking asymptotic statistical properties 
corresponding to taking the limit T — > oo. Let : 



7Ta(.) = Hm 7T,- 



(T) 



(■), 



(2) 



be the empirical measure for the raster plot Q, then 
v = ttq for z/-almost every u>. 

Gibbs measures. The explicit form of v is not 
known in general. Moreover, one has experimen- 
tally only access to finite time raster plots and the 
convergence to tt^ in (fJJ can be quite slow. Thus, 
it is necessary to provide some a priori form for the 
probability v. We call a statistical model an ergodic 
probability measure v' which can serve as a proto- 
typal for v. Since there are many ergodic measure, 
there is not a unique choice for the statistical model. 
Note also that it depends on the parameters I\ We 
want to define a procedure allowing one to select a 
statistical model. 

A Gibbs measure with potential tp 
is an ergodic measure v^. r such that 



y i/>;r 



0" 



)1,T 



where Z^\ip) = 



da 



with respect to v^,. r . 



(3) 



the average value of <fi with respect to v^- <r . The 
Gibbs measure obeys a variational principle. Let v 
be a <7 r -invariant measure. Call h(v) the entropy of 
v. Let Ai™ v be the set of invariant measures for a r , 
then: 



sup 



(h(u) + v(ip)) = h(^. r )+^. r (ip). 

(4) 



Adaptation rules. We consider adaptation mecha- 
nisms where synaptic weights evolve in time accord- 
ing to the spikes emitted by the pre- and post- synap- 
tic neuron. The variation of Wij at time t is a func- 
tion of the spiking sequences of neurons i and j from 
time t — T to time t, where T is time scale charac- 
terizing the width of the spike trains influencing the 
synaptic change. In this paper we investigate the ef- 
fects of synaptic plasticity rules where the character- 
istic time scale T is quite a bit larger than the time 
scale of evolution of the neurons. 

On practical/numerical grounds we proceed as 
follows. Let the neurons evolve over a time win- 
dows of width T, called an "adaptation epoch" dur- 
ing which synaptic weights are constant, and record 
the spikes trains. From this record, we update the 
synaptic weights and a new adaptation epoch begins. 
We denote by t the update index of neuron states in- 
side an adaptation epoch, while r indicates the up- 
date index of synaptic weights. Neuronal time is re- 
set at each new adaptation epoch. 

As an example we consider here an adaptation 
rule of type: 



(r) 



y 
(A 



(r+l) 



(t) 



l ) W ij + T=2c9ij(N 0iT , Mo, T ) 



The first term parametrized by A 6 
sive LTD while: 



[0,1" 



" (5) 
mimics pas- 



J^„(T) exp J2t=i t l J ( a r 1 ^)' ^r- 7 ^ being the set of ad- 
r 

missible spiking patterns sequences of length T. 
Note that for the forthcoming discussion it is im- 
portant to make explicit the dependence of these 
quantities in the parameters T. 

The topological pressur^is the limit P r (ip) = 
lirriT^oo h log(-Z'| T ' ) If 4> is another potential 

2 e.g. by minimising the Kullback-Leiber divergence between the empirical measure and v 1 . 

3 Fix < < 1. Define a metric on S r by (cD, u>') = © p , where p is the largest integer such that u>(t) = <jj(t'), 1 < t < p — 1. For 
a continuous function i/> : S r — » R and n > 1 define var n ip = sup {|i/>((2>) — ip(uj')\ : u>(t) = <*>'(£), 1 < i < n}. ip is called a potential 
if var n {ij)) < C@ n , n = 1, 2 . . . , where C is some positive constant 1* 351 . 

4 This quantity is analogous to thermodynamic potentials in statistical mechanics like free energy or the grand potential Q = U — TS — 
fiN = —PV, where P is the thermodynamic pressure. Equation {5) expresses that Pr(ip) is the generating function for the cumulants of the 
probability distribution. Equation (4) relates Gibbs measures to equilibrium states (entropy maximisation under constraints). 



J0,T • 



No,r) = 



T-C 



C 



t=C u=-C 



f(u)u>i(t+u)u>j(t), 



with : 



/(*) = 



(6) 



A_e? 

A + e 

0, 



x < 0; 
x > 0; 
a; = 0, 



and with A- < 0, A + > 0, provides an example 
of STDP implementation. Note that since f(x) be- 
comes negligible as soon as x S> r + or x <C — r_, 
we may consider that f(u) = whenever \u\ > 

C = f 2 max(r + , r_). The parameter e > is cho- 
sen small enough to ensure adiabatic changes while 
updating synaptic weights. 

Updating the synaptic weights has several con- 
sequences. It obviously modifies the structure of the 
synaptic graph with an impact on neurons dynamics. 
But, it can also modify the statistics of spikes trains 
and the structure of the set of admissible raster plots. 
Let us now discuss these effects within more details. 

Dynamical effects of synaptic plasticity. Set 

S i{&) = Eu=-c/( M W") and H ij(v) = 
ujj(0)Si(u>). Then the adaptation rule writes, as 

T — > oo, in terms of the empirical measure 7T S ( T ) : 

SWg) = e ((A - 1)WV + W [Hij]) , (7) 

where Cj^ is the raster plot produced by neurons in 
the adaptation epoch r. 

The term Si(u>) can be either negative, induc- 
ing LTD or positive inducing LTP. In particular, its 
average with respect to the empirical measure ttq( T ) 
reads: 

7r aW (S i ) = /9r i (r) (8) 

where: 

r c c -i 



1 - e T - 1 - e T + 

(9) 

and where r<(r) = w^w (u>i) is the frequency rate of 
neuron i in the r-th adaptation epoch. 

The term /3 neither depend on uj nor on 
r, but only on the adaptation rule parameters 
A-, A + , r_,r+, C. Equation $Q leads 3 regimes. 

• Cooperative regime. If /3 > then 

^tDt-r) (<%) > 0. Then synaptic weights have 
a tendency to become more positive. This cor- 
responds to a cooperative system |36| . In this 
case, dynamics become trivial with neurons 
firing at each time step or remaining quiescent 
forever. 

• Competitive regime. On the opposite if /? < 
synaptic weights become negative. This cor- 
responds to a competitive system (36]. In this 
case, neurons remain quiescent forever. 

• Intermediate regime. The intermediate 
regime corresponds to p ~ 0. Here no clear 
cut tendency can be distinguished from the av- 
erage value of Si and spikes correlations have 
to be considered as well. In this situation, the 



effects of synaptic plasticity on neurons dy- 
namics strongly depend on the values of the 
parameter A. A detailed discussion of these ef- 
fects will be published elsewhere. We simply 
provide an example in Fig. Q] 

Variational form. We now investigate the effects of 
the adaptation rule (O on the statistical distribution 
of the raster plot. To this purpose we use, as a sta- 
tistical model for the spike train distribution in the 
adaptation epoch t, a Gibbs distribution i/^m with 
a potential ip^ and a topological pressure 
Thus, the adaptation dynamics results in a sequence 
of empirical measures {tt^m } t>0 , and correspond- 
ing statistical models {^m} >0 > con "esponding to 
changes in the statistical properties of raster plots. 

These changes can be smooth, in which case 
the average value of observables changes smoothly. 
Then, using (O we may write the adaptation rule © 
in the form: 

= M 

= V W=W(T ,P + (W- W«).0(W (T) . •)] ' 

(10) 

where is the matrix of components <j>ij with 0y = 
4>ij (Wij , u>) = e ((A — l)Wij + Hij ) and where we 
use the notation (W - W^).^ = £f,=i(Wy - 

This amounts to slightly perturb the potential 
t/i (t) with a perturbation (W - W (r) ).0. In this 
case, for sufficiently small e in the adaptation rule, 
the variation of pressure between epoch r and r + 1 
is given by: 

5P {T) - <5W (T) .^(r) [<f>] = 5W {T) .5W {T) > 

Therefore the adaptation rule is a gradient system 
which tends to maximize the topological pressure. 

The effects of modifying synaptic weights can 
also be sharp (corresponding typically to bifurca- 
tions in (Q])). We observe in fact phases of regular 
variations interrupted by singular transitions (see fig. 
[TJi. These transitions have an interesting interpreta- 
tion. When they happen the set of admissible raster 
plots is typically suddenly modified by the adapta- 
tion rule. Thus, the set of admissible raster plots ob- 
tained after adaptation belongs to E r ( T ) (~l S r ( T +i). 
In this way, adaptation plays the role of a selective 
mechanism where the set of admissible raster plots, 
viewed as a neural code, is gradually reducing, pro- 
ducing after n steps of adaptation a set n" =1 £ r M 
which can be small (but not empty). 

If we consider the situation where (|T|i is a neu- 
ral network submitted to some stimulus, where a 
raster plot u> encodes the spike response to the stim- 
ulus when the system was prepared in the initial 
condition X — * ui, then E r is the set of all possible 



raster plots encoding this stimulus. This illustrates 
how daptation results in a reduction of the possible 
coding, thus reducing the variability in the possible 
responses. This property has also been observed in 

EH [32 HH [301 

Equilibrium distribution. An "equilibrium point" 
is defined by SW — and it corresponds to a maxi- 
mum of the pressure with : 

\4>] = o, 

where v* is a Gibbs measure with potential ip*. This 
implies: 

w* = —L-z,* [ft] (ii) 

which gives, component by component, and making 
H explicit: 



w *j = 73T £ f{u)v * ' (12) 

u=-C 

where v* [lu 3 (O)wi(it)] is the probability, in the 
asymptotic regime, that neuron i fires at time u and 
neuron j fires at time 0. Hence, in this case, the 
synaptic weights are purely expressed in terms of 
neurons pairs correlations. As a further consequence 
the equilibrium is described by a Gibbs distribution 
with potential <p. 

Therefore, if the adaptation rule (0 converges^ 
to an equilibrium point as r — > oo, it leads naturally 
the system to a dynamics where raster plots are dis- 
tributed according to Gibbs distribution, with a pair 
potential having some analogy with an Ising model 
Hamiltonian. This makes a very interesting link 
with the work of Schneidman et al [37] where they 
show, using theoretical arguments as well as em- 
pirical evidences on experimental data for the sala- 
mander retina, that spike trains statistics is likely de- 
scribed by a Gibbs distribution with an "Ising" like 
potential (see also |38|). 

3 Conclusion. 

In this paper, we have introduced a mathematical 
framework where spikes trains statistics, produced 
by neural networks, possibly evolving under synap- 
tic plasticity, can be analysed. It is argued that Gibbs 
distribution, arising naturally in ergodic theory, are 
good candidates to provide efficient statistical mod- 
els for raster plots distribution. We have also shown 
that some plasticity rules can be associated to a vari- 
ational principle. Though only one example was pro- 
posed, it is easy to extend this result to more general 

5 This convergence is ensured for sufficiently small e. 



examples of adaptation rules. This will be the subject 
of an extended forthcoming paper. 

Concerning neuroscience there are several ex- 
pected outcomes. The approach proposed here is 
similar in spirit to the variational approaches dis- 
cussed in the introduction (especially 1221 l23l [24]). 
Actually, the quantity L minimized in ll23l l24l can 
be obtained, in the present setting, with a suitable 
choice of potential. But the present formalism con- 
cerns neural networks with intrinsic dynamics (in- 
stead of isolated neurons submitted to uncorellated 
Poisson spikes trains). Also, more general models 
could be taken into account (for example Integrate 
and Fire models with adaptive conductances |34|). 

A second issue concerns spike trains statistics. 
As claimed by many authors and nicely proved by 
Schneidman and collaborators 1371 l38l . more elab- 
orated statistical models than uncorrected Poisson 
distributions have to be considered to analyse spike 
trains statistics, especially taking into account cor- 
relations between spikes [39|, but also higher order 
cumulants. The present work show that Gibbs dis- 
tributions, obtained from statistical inference in real 
data by Schneidman and collaborators, may naturally 
arise from synaptic adaptation mechanisms. We also 
obtain an explicit form for the probability distribu- 
tion depending e.g. on physical parameters such as 
the time constants t + , t_ or the LTD/LTP strength 
A + , A- appearing in the STDP rule. This can be 
compared to real data. Also, having a "good" statis- 
tical model is a first step to be able to "read the neural 
code" in the spirit of [17], namely infer the condi- 
tional probability that a stimulus has been applied 
given the observed spike train, knowing the con- 
ditional probability that one observes a spike train 

given the stimulus. 

Finally, as a last outcome, this approach opens 
up the possibility of obtained a specific spike train 
statistics from a deterministic neurons evolution with 
a suitable synaptic plasticity rule (constrained e.g. by 
the potential 0). 
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Figure 1 . Effect of synaptic plasticity on: (up left) synaptic weights histogram; (up right) Period of the orbit; (down 
left) some synaptic weights evolution; (down right) raster plot of a few neurons. 



